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SUMMARY/OVERVIEW: 

The  objective  of  this  study  is  the  fundamental  understanding  of  fuel  disintegration  and  mixing  in 
a supercritical  environment  (relative  to  the  fuel)  in  order  to  determine  parameter  regimes 
advantageous  to  mixing.  The  approach  is  based  on  developing  a model  of  a supercritical, 
turbulent  jet  mixing  with  surrounding  fluid.  The  method  is  one  that  combines  the  modeling  of 
supercritical  fluids  with  a systematic  development  based  on  the  Large  Eddy  Simulation  (LES) 
approach.  This  systematic  development  includes  a consistent  protocol  based  upon  Direct 
Numerical  Simulations  (DNS)  for  developing  a Subgrid  Scale  Model  (SGS)  appropriate  to 
supercritical  fluids,  rather  than  choosing  in  an  ad  hoc  manner  an  existing  SGS  model  developed 
under  assumptions  inconsistent  with  supercritical  fluid  behavior.  This  SGS  model  will  be  used  in 
the  LES  of  a supercritical  turbulent  jet. 

TECHNICAL  DISCUSSION 

The  first  step  in  this  effort  consisted  in  the  development  of  a DNS  and  associated  code  for 
a supercritical  shear  layer.  The  shear  layer  is  chosen  as  the  simplest  configuration  from  which 
one  may  obtain  results  that  will  establish,  through  DNS,  a framework  for  understanding  mixing 
of  supercritical  fluids  at  small  scales.  Because  the  turbulent  behavior  at  the  Kolmogorov  scales  is 
universal  (i.e.  geometric-configuration  independent),  once  this  behavior  is  determined  (e.g.  in  a 
shear  layer  geometric  configuration),  it  can  be  utilized  to  model  a variety  of  geometric 
configurations  (for  the  same  fluids)  having  the  same  physics.  Since  the  interest  of  AFOSR  is  in 
hydrocarbons,  the  chosen  shear  layer  is  that  of  heptane  (lower  stream)  in  nitrogen.  Having 
conducted  preliminary  investigations  in  which  the  layer  did  not  reach  transition  (Miller  et  al., 
2001),  and  having  analyzed  the  causes  for  the  lack  of  transition  (Okong’o  and  Bellan,  2000a),  the 
activities  were  further  focused  on  obtaining  one  or  several  transitional  states  and  analyzing  their 
behavior.  To  this  end,  there  were  several  areas  of  emphasis: 

(i)  The  real-gas  boundary  conditions  of  Okong’o  and  Bellan,  2000b  (OBOOjcp)  were 
utilized  in  the  DNS  code  and  two  simulations  were  conducted  differing  only  by  the  value  of  the 
initial  Reynolds  number;  the  wavelength  of  the  initial  forcing  was  the  typical  one  from  the 
incompressible  stability  analysis.  The  results  were  documented  in  Okong’o  and  Bellan,  2000c 
(OBOOjfrn). 
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The  transitional  state  was  identified  by  global  manifestations  such  as  rapid  and  sustained 
momentum  thickness  growth,  high  rate  of  positive  spanwise  vorticity  evolution,  increased 
enstrophy,  large  product  thickness,  and  large  momentum  thickness  based  Reynolds  number. 
Analysis  of  the  vorticity  and  vorticity  magnitude  budgets  at  transition  highlighted  the  dominance 
of  the  stretching  and  tilting  effect  in  the  production  of  spanwise  vorticity,  and  of  the  viscous 
contribution  in  the  vorticity  magnitude,  although  at  some  crosstream  locations  stretching  and 
tilting  exceeded  viscous  effects. 

Visualizations  of  the  spanwise  vorticity  at  the  transitional  state  revealed  considerable 
local  positive  vorticity  both  in  the  braid  and  between  the  braid  planes.  Furthermore, 
visualizations  of  the  magnitude  of  the  density  gradient  displayed  concentrated  regions  of  high 
values  exhibiting  a convoluted  and  distorted  appearance;  these  regions  may  be  optically 
identifiable  owing  to  the  change  in  the  index  of  refraction.  Given  the  commonality  of  distinctive 
optical  features  between  temporal  and  spatial  shear  layers,  the  observation  of  such  regions  in 
supercritical  jet  experiments  provided  encouragement  that  the  essence  of  the  situation  has  been 
captured.  Analysis  of  the  fluid  composition  in  these  regions  of  high  density  gradient  magnitude 
revealed  that  the  predominant  species  is  that  in  the  lower,  entrained  stream  with  small  amounts  of 
the  entraining  fluid  transported  and  dissolved  into  it.  Moreover,  considerations  based  on  the 
value  of  the  mass  diffusion  factor  identified  these  regions  to  contain  a highly  non-ideal  mixture, 
implying  mixing  difficulties  at  the  molecular  level  resulting  from  the  thermodynamic  properties 
of  the  mixture.  Visualizations  of  the  compression  factor  displayed  large  departures  from  perfect 
gas  behavior  in  the  lower  heptane  stream  as  well  as  in  the  layer. 

The  issue  of  the  thermodynamic  state  of  the  mixture  at  transition  was  also  explored 
through  examination  of  departures  from  the  critical  locus  both  in  the  braid  and  between  braid 
planes.  The  fluid  was  supercritical  everywhere,  as  the  temperature  was  always  larger  than  that  of 
the  local  critical  point.  However,  the  pressure  was  above  the  local  critical  point  outside  and 
below  the  local  critical  point  inside  the  vortices,  with  the  critical  locus  mapping  the  coherent 
vortices.  The  considerable  reduction  in  the  difference  between  the  local  values  of  the  pressures 
and  critical  pressure  inside  the  vortices  was  due  to  the  change  in  fluid  composition. 

To  investigate  the  primary  mechanisms  responsible  for  dissipation  both  during  the 
evolution  of  the  layer  and  at  transition,  the  irreversible  entropy  production  was  calculated  and 
analyzed.  The  three  contributions  to  the  dissipation  arising  from  viscous  stresses,  species  fluxes 
and  heat  fluxes  were  calculated  and  compared.  Volume  averages  as  well  as  RMS  were  computed 
at  time  stations  after  the  second  pairing,  but  prior  to  transition,  and  at  transition.  Moreover,  by 
filtering  the  DNS  solution  and  calculating  the  difference  between  the  filtered  and  unfiltered 
variables,  the  contribution  of  the  small  scales  (SGS)  was  evaluated.  The  overwhelming  SGS 
contribution  both  to  the  average  and  the  RMS  is  due  to  viscous  effects,  and  the  SGS  activity  is 
also  dominated  by  viscous  effects. 

Since  the  transitional  state  was  dominated  by  viscous  effects,  one  could  argue  that  the 
above  conclusions  are  universal,  i.e.  independent  of  the  species  and  only  governed  by  the 
supercritical  state.  Such  an  argument  may  be  in  error  since  some  of  the  important  aspects  of  the 
transitional  state,  such  as  the  fluid  composition  in  the  high  magnitude  density  regions  may  be  a 
result  of  the  mixture  thermodynamics.  Ascertaining  the  species  independent  aspects  of 
supercritical  mixing  layer  behavior  will  constitute  one  of  our  future  priorities. 

Assumed  PDF  modeling,  of  particular  interest  for  reacting  flow  representations,  was 
shown  to  have  only  limited  potential.  Neither  the  beta-density  nor  the  Gaussian  PDF,  both 
constructed  using  the  DNS  calculated  moments,  predicted  reliably  the  mixing  process  or  the 
transitional  state,  respectively.  The  same  result  was  found  at  the  subgrid  scale.  Moreover,  the 
temperature  and  partial  density  PDFs  were  well  (negatively)  correlated,  indicating  that  the 
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possible  modeling  of  the  reaction  rate,  which  is  a joint  temperature-species  PDF,  by  the  product 
of  the  marginal  PDFs  holds  poor  prospects.  This  good  correlation  was  obtained  both  for  the  DNS 
scale  PDF  and  the  SGS  scale  PDF  (the  filtered  density  function  of  Pope,  1991).  Finally,  it  should 
be  mentioned  that  the  conclusions  regarding  the  PDF  representation  are  the  result  of  studying 
only  two  realizations.  It  was  realized  that  additional  work  is  necessary  to  obtain  more  definitive 
conclusions. 

(ii)  The  investigation  described  in  (i)  indicated  the  need  for  additional  realizations  of 
transitional  states.  However,  each  DNS  simulation  is  very  computationally  intensive,  and  thus  a 
good  strategy  was  deemed  necessary  to  minimize  computational  costs.  Therefore,  a linear 
inviscid  stability  analysis  for  real  gas  was  conducted  with  the  goal  of  finding  other  forcing 
wavelengths  that  perhaps  promote  quicker  transition.  The  results  of  the  stability  analysis  and  the 
ensuing  simulations  were  documented  in  Okong’o  and  Bellan,  2001.  Two-dimensional  (2D) 
stability  results  showed  that  the  unstable  growth  rates  for  the  compressible  flow  are  smaller  at  all 
wavenumbers  than  the  equivalent  growth  rates  found  from  an  incompressible  analysis;  also,  the 
most  unstable  wavelengths  occurred  at  smaller  wavenumbers  with  increasing  density 
stratification.  2D  mixing  layer  simulations  showed  that  at  the  same  initial  mean  flow  conditions 
the  solution  is  virtually  independent  of  the  excitation  wavelength.  This  observation  raised  the 
alluring  prospect  of  finding  the  shortest  streamwise  compressible  unstable  wavelength  at  which 
rollup  and  pairing  could  be  obtained  in  three-dimensional  (3D)  simulations;  for  such  a 
wavelength,  the  size  of  the  domain  could  be  reduced  (since  it  is  four  times  the  perturbation 
wavelength)  leading  to  a reduced  number  of  grid  points  for  the  same  resolution.  However,  in 
contrast  with  the  2D  problem  where  for  each  wavelength  there  is  a single  eigenvalue  (i.e. 
propagation  speed),  in  3D  it  was  found  that  additional  to  the  wavelength  of  the  perturbation,  the 
angle  of  the  wave  direction  in  the  streamwise-spanwise  plane  must  be  specified  as  well.  To 
minimize  the  wavelength  of  the  spanwise  perturbation,  it  was  shown  that  one  must  maximize  the 
angle  representing  the  direction  of  the  perturbation.  For  the  purpose  of  applications  to  3D  flow 
simulations,  this  implies  that  shorter  wavelengths  in  the  streamwise  direction  are  counter 
balanced  by  the  need  to  use  a longer  wavelength  in  the  spanwise  direction.  Therefore,  to  avoid 
enlarging  the  spanwise  domain  size,  the  heuristic  approach  of  using  non-eigenfunction 
perturbations  was  adopted  and  simulations  conducted  for  five  shear  layer  were  compared;  for 
each  of  these  simulations  the  streamwise  domain  size  was  four  times  the  perturbation  wavelength 
and  the  calculation  was  initiated  with  four  vortices  leading  after  two  pairings  to  an  ultimate 
vortex.  Among  these  five  simulations,  two  layers  (differing  by  the  value  of  the  initial  Reynolds 
number)  were  excited  with  the  incompressible  most  unstable  wavelength.  The  third  and  fourth 
layers  (differing  by  the  value  of  the  initial  Reynolds  number)  were  excited  with  the  shortest 
compressible  unstable  wavelength  from  the  2D  stability  analysis  and  the  same  amplitudes  as  the 
first  two  layers.  The  fifth  layer  was  excited  with  the  shortest  compressible  unstable  wavelength 
from  the  two-dimensional  stability  analysis  and  a larger  amplitude  of  the  3D  excitation.  Only 
three  of  these  five  layers  achieved  transition,  and  the  expectation  of  reduced  CPU  time  for  the 
layer  perturbed  with  the  shortest  compressible  unstable  wavelength  did  not  materialize.  Global 
evolution  and  contour  plots  showed  that  the  inability  of  transition  achievement  by  two  of  the 
layers  perturbed  with  the  shortest  compressible  unstable  wavelength  is  due  either  to  the  inability 
of  the  formed  small  scales  to  produce  substantial  vorticity  or  to  the  early  formation  of  small 
scales  and  vorticity  production  which  destroy  the  coherence  of  the  ultimate  vortex  preventing  its 
growth  and  initiating  its  decay. 

The  results  showed  that  transitional  states  with  the  same  momentum  thickness  based 
Reynolds  number  can  be  achieved  by  perturbing  with  a variety  of  unstable  wavelengths,  and  that 
they  require  similar  CPU  time.  Simulations  performed  with  a shorter  wavelength  perturbation 
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require  larger  initial  Reynolds  numbers  to  reach  a transitional  state.  Transitional  states  obtained 
with  different  perturbation  wavelengths  had  different  global  and  detailed  characteristics,  and 
therefore  such  simulations  were  recommended  for  enlarging  the  database  necessary  to  extract 
turbulence  models  utilizable  in  the  mathematical  description  of  turbulent  spatial  flows. 

Further  investigations  are  now  devoted  to  enlarging  the  transitional  states  database  by  (1) 
considering  other  binary  system  species,  particularly  those  for  which  contemporary  experimental 
information  is  available  (e.g.  LOX/He),  and  (2)  considering  additional  perturbation  wavelengths 
to  obtain  transitional  states  with  different  characteristics. 
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Figure  1 Density  Gradient  Magnitude  at  t*=135  (corresponding  to  the  attainment  of  a transitional  state),  in 
the  braid  plane  (left)  and  in  the  between-the-braid  plane  (right).  The  initial  Reynolds  number  is  600,  the 
initial  Mach  number  is  0.4,  and  the  freestream  temperatures  are  1000  K for  the  upper,  nitrogen,  stream  and 
600  K for  the  lower,  heptane,  stream. 
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